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ABSTRACT 


The  determination  of  periodic  edge  waves  along  a  straight  infinite  beach 
can  be  formulated  as  a  2-dimensional  eigenvalue  problem  for  a  generalized 
wedge.  A  least  squares  boundary  collocation  method  for  computing  the  lower 
discrete  eigenvalues  is  described,  with  an  error  estimate  for  the  least 
eigenvalue.  Numerical  examples  are  presented. 


AMS(MOS)  Subject  Classification:  65N25,  65N35 

Key  Words:  Boundary  collocation.  Edge  waves.  Eigenvalue  problem.  Sloping 
beach 

Work  Unit  Number  3  -  Numerical  Analysis  and  Scientific  Confuting 
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SIGNIFICANCE  AND  EXPLANATION 


Computational  eigenvalue  problems  for  unbounded  domains  often  present 
difficulties  in  execution.  Thus  it  is  of  interest  to  demonstrate  the 
viability  of  the  boundary  collocation  method  for  problems  on  unbounded  domains 
for  which  few  analytic  results  are  known. 

For  demonstration  the  Stokes-Ursell  edge  waves  for  general  bottom 
profiles  are  examined.  Even  though  other  methods  are  available  for 
computation  of  these  eigenvalues,  the  flexibility  of  the  present  method  might 
be  an  advantage  in  technical  problems.  It  should  also  be  applicable  to 
different  problems  on  semi-infinite  domains. 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC,  and  not  with  the  authors  of  this  report. 


EDGE  WAVES  BY  BOUNDARY  COLLOCATION 


L.  Reichel  and  P.  Chapman 


1.  Introduction 

^  The  aim  of  this  paper  is  to  describe  a  computational  method  for  deter¬ 
mining  traveling  edge  waves  along  a  straight  infinite  beach  with  a  cylin¬ 
drical,  not  necessarily  plane  bottom. '^Assuming  the  waves  are  periodic  along 
the  shore  line  (the  z-axis),  this  problem  can  be  reduced  to  a  2-dimensional 
eigenvalue  problem  on  a  generalized  wedge,  V,  in  the  x,y-plane,  see  Ursell 
[9].  The  top  boundary  of  the  wedge,  representing  the  (horizontal)  free 
surface  of  the  fluid,  is  the  positive  x-axis.  The  other  boundary  of  the  wedge 
is  below  the  positive  x-axis  and  represents  the  sea  bottom,  see  Figure  1.1. 

The  y-coordinate  of  the  sea  bottom,  y  *  b(x),  is  assumed  to  be  a  smooth 
(C2)  function  of  x,  with  b(0)  *  0,  b(x)  <  0,  x  >  0  and  with  a  finite, 

negative  Lim  -  . 

x+«  X 
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Figure  1 . 1 
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The  traveling  waves  of  interest,  see  Ursell  [9],  p.  83,  are 

1>(x,y,z,t)  :■  cos(ks)^(x,y)eipt  ,  (1.1) 

where  $  and  the  real  constant  p  satisfy 

A  ♦  -  ♦  on  V,  (1.2) 

xy 

+  =*:  w2^  on  top:  y=0,  x>0,  (1.3) 

an  Kg 

^  m  0  on  bottom:  y=b(x),  x>0,  (1.4) 

an 

E  *-  ii(lX*|2  +  ^2)  «3S  <  -  .  (1.5) 

V 

n  denotes  the  outward  unit  normal,  g  is  the  gravitational  acceleration 

and  k  is  the  same  real  constant  as  in  (1.1).  The  dimensionless  problem 

eigenvalue  is  denoted  by  w  .  Ursell  [9]  studied  the  case  of  a  straight 

bottom  (i.e.  b(x)  ■  cx  )  and  found  a  finite,  simple,  discrete,  spectrum  with 

eigenfunctions  $  of  finite  energy  E  (equation  1.5),  and  a  continuous 

2 

spectrum  (of  non-finite  energy)  above  a  cut  off  value  of  w  . 

Ursell  [9]  further  observed  that  the  number  of  eigenpairs  of  the  problem 
(1.2-1. 5)  increases  as  the  interior  angle  of  the  wedge  decreases.  A  recent 
analysis  based  on  shallow  water  theory  of  Chapman  [1]  suggests  that  this  angle 
is  also  a  controlling  parameter  when  the  bottom  is  not  straight. 

In  the  present  work  the  central  interests  are  the  following 
1.  To  study  the  use  of  a  Boundary  Collocation  Method  (BCM)  in  the 
computation  of  eigenpairs. 

To  investigate  the  role  of  the  wedge  angle  as  a  dominant  parameter. 


2 


The  attraction  of  the  BCM  lies  in  its  applicability  to  a  variety  of  domains, 
its  ease  of  programming,  and  its  conceptual  simplicity.  The  method  reduces 
the  compute t iona 1  problem  to  a  1-dimensional  problem.  Murashima  et  al.  [4] 
call  the  BCM  the  charge  simulation  method  (CSM)  and  state  "The  charge 
simulation  method  [...]  is  often  used  in  some  part  of  Electrical  Engineering 
as  an  efficient  way  for  solving  Laplace's  equation.  The  CSM,  is  not 
considered  to  be  straight  forward  because  contour  points  t=  collocation 
points]  and  source  points  I*  singular  points  of  the  basis  functions]  must  be 
determined  by  trial  and  error."  We  eliminate  the  guessing  in  the  CSM  by 
presenting  a  strategy  for  selecting  collocation  points  when  the  set  of  basis 
functions  is  given.  We  will  choose  the  number  of  collocation  points  so  that 
an  overdetermined  system  of  equation  is  obtained.  This  we  solve  by  a  Galerkin 
method,  i.e.  the  boundary  collocation  method  to  be  described  is  a  Galerkin- 
coliocation  method.  For  each  set  of  basis  functions  the  determination  of 
eigenpairs  is  quite  rapid,  and  we  carry  out  the  confutation  for  several 
systematically  chosen  sets  of  basis  functions.  A  simply  computable  error 
estimate  for  the  lowest  eigenvalue  reveals  which  set  of  basis  functions  yields 
best  accuracy.  This  makes  it  unnecessary  for  the  user  to  guess  good  source 
points.  The  idea  to  iteratively  choose  source  points  has  previously  been 
advocated  by  Mathon  and  Johnson  [2]  for  Dirchlet  problems.  The  number  of 
distinct  source  points  we  need  is,  however,  much  smaller  than  in  the  examples 
in  [2],  This  is  important  since  the  computational  effort  increases  rapidly 
with  the  number  of  distinct  source  points,  see  section  2. 

This  paper  is  organized  as  follows:  Section  2  contains  the  description 
of  the  BCM,  Section  3  contains  numerical  examples,  and  in  Section  4  we  derive 
an  easily  computable  error  estimate  for  the  least  eigenvalue.  The 
eigenproblem  under  consideration  has  previously  been  studied  by  methods  of 
geometrical  optics,  see  Meyer  [3]  and  Shen  et  al.  [8]. 


rVi_  » 


.  A  Boundary  Collocation  Method 


He  first  describe  the  allocation  of  collocation  points  for  a  given  set  of 
eigenfunctions.  This  description  is  slightly  more  general  than  needed  for  the 
numerical  examples  in  sections  3,  because  we  want  generalizations  to  other 
problems  to  be  simple.  Presently,  we  only  can  motivate  the  allocation 
heuristically,  but  numerical  examples  show  that  the  allocation  is  appropriate. 

Let  {(Xj,  y  )}  ^  be  a  point  set  exterior  to  V  U  9P .  For  an  arbitrary 
point  (x,y),  let  (r^,  6^)  be  the  polar  coordinates  for  (x-x^,  y-y^), 
j  *  1(1)£.  We  will  approximate  the  eigenfunctions  of  (1.1)— (1.5)  by  linear 
combinations  of  the  functions 


Kfl(r^),  Kj  (r^  )cos(6j  ),  (r^  )sin(  0^ ) ,  K2  (r  _.)cos(  20^  ) , . . . 


,,K  (r  )cos(pfl . ),  K  (r  )sin(p6  ),  j=1(1)i, 

pj  3  Pj  j 


(2.1) 


for  some  integer  p  >  0.  K^(r)  denotes  the  j  order  modified  Bessel 

function  of  the  2nd  kind.  Any  linear  combination  of  functions  (2.1)  satisfies 
(1.2)  and  (1.5).  Conditions  (1.3)  and  (1.4),  we  impose  at  collocation  points 


We  turn  to  the  allocation  of  collocation  points.  We  are  led  to  the 


allocation  strategy  by  regarding  approximation  of  analytic  functions  f(w)  on 


regions  ft  in  the  extended  complex  plane  C.  The  complement  of  ft,  ft  ,  is 

c 

assumed  to  be  open,  bounded,  simply  connected,  and  have  a  smooth  boundary 

l 

3fl.  Let  {Wj}_j=1  be  a  finite  set  of  distinct  points  in  Q^,  and 
approximate  f(w)  by  a  linear  combination  of  the  functions 


1,  (w-w ^ \  (w-Wj)  2,...,  (w-Wj )  p,  j=1(1)£.  (2.2) 


•  •  ••  .vy, 


The  linear  combination  is  to  be  determined  by  interpolation  at  points  on  9ft. 
Regard  w^,  unit  point  charges,  which  generate  a  potential 

“7  1  *n|w-w.  |.  Let  o(w)  be  the  charge  distribution  on  3ft  of  total  charge 

Vi  j 

1,  so  that  o(w)  makes  3ft  an  equipotential  curve  in  the  presence  of  the 

A 

point  charges  at  Wj.  It  is  known  that  m  =  pi+1  interpolation  points 

(■  collocation  points)  on  for  determining  a  linear  combination 

c 

of  functicrs  (2.2)  should  be  equidistributed  w.r.t.  c(w),  i.e.  w1  is 
c 

arbitrary  and  w  ^  ,  j  >  1,  should  satisfy 


w 

/k+1a(w)|dw|  ■  — ,  k*1(1)m-1. 


(2.3) 


where  integration  is  done  along  3ft  in  postive  direction.  For  details,  see 
Reichel  [6].  o(w)  can  be  determined  by  solving 


q+  /  in|w-w|o(w) |dw|  =  7  I  tn|w-w  |,  w  e  3ft 
9ft  *j=1  5 


(2.4) 


/  7(w)|dw|  =  1, 

9ft 


where  q  is  a  constant  also  to  be  determined,  see  [6]  or  [7J.  Let  F(w)  be 
a  conformal  map  from  ft^  to  the  complement  of  V  U  3P.  Let 


x.  +  iy .  :=  F(w . ) ,  j=1(1)fc. 
3  3  3 


F  can  be  continued  to  a  bijective  map  on  ft  U  9 ft,  and  we  let  the  collocation 
c  c 

points  (x_j,  y_.),  j=1(1)m  be  defined  by 


x^  +  iy°  =  F(w^),  j  =  1(1 )m. 


(2.5) 


-L*  ■.*  «.-V  V*  •;  •.*  %•  s’  •  •.  . 


We  motivate  this  choice  by  the  invariance  of  the  potential  under  conformal 


mapping,  and  by  the  fact  that  for  some  constants  c  *  0 , 

P 

V  V 


Kg  (r j  )  +  Cgin(rj)  as  r..  ♦  0, 


independent  of 


-ip6.  -pe-ip0. 

K  (r,)e  J  ♦  c  r.  J  as  r.  0,  p  =  1,2,3,..., 

P  3  P  3  3 


which  shows  that  the  functions  (2.1)  behave  essentially  like  harmonic 
functions  as  r^  ♦  0,  j=1(1)i.  For  harmonic  functions  the  same  allocation  as 
for  analytic  functions  is  suitable. 

The  application  to  problem  (1.1)- (1.5)  permits  some  simplifications. 

c  c 

When  allocating  the  collocation  points  (x_.,  y  ),  we  approximate  the 
generalized  wedge  V  by  a  wedge  D  with  straight  sides  and  with  vertex  at 
the  origin.  Let  the  interior  angle  of  D  at  the  origin  be  0x.  Turn  D  to 
be  symmetric  w.r.t.  the  x-axis  and  so  that  D  has  the  negative  x-axis  in  its 
interior.  A  conformal  mapping  F  which  maps  the  exterior  of  the  unit  circle 
on  the  complement  of  D  is  explicitly  known 


(2.6) 


The  charge  density  function  o(w(s))  we  seek  in  the  form  of  a  trigonometric 

A 

polynomial  of  degree  n.  When  3£)  is  a  circle  then  this  approximation  of 

A  A 

o(w(s))  can  be  determined  in  o(n  in(n))  operations,  see  [5]  especially 
lemma  2  for  details.  We  choose  the  angle  of  D  so  that  properly  turned,  D 
approximates  V  well  near  the  origin.  Having  determined  the  collocation 
points  on  3d,  we  map  them  onto  3P  by  some  mapping  which  approximately 
preserves  distance  from  the  vertex. 


In  the  computations  we  have  chosen  the  number  collocation  points,  m,  as 


2-3  times  n,  the  number  of  basis  functions  (2.1),  i.e.  n  =  (2p  +  1)4.  The 
allocation  is  good  in  the  sense,  that  doubling  the  number  of  collocation 
points  has  in  no  case  when  checked  changed  the  order  of  magnitude  of  the  error 
in  the  eigenvalues. 

We  now  turn  to  the  computation  of  eigenvalues  and  eigenfunctions  with  a 

given  set  of  basis  functions  (2.1).  First  we  allocate  m  collocation  points 

as  described,  and  make  sure  one  point  is  placed  at  the  vertex.  The  vertex 

point  is  regarded  as  2  points,  since  we  will  impose  both  conditions  (1.2)  and 

(1.3)  at  the  vertex.  Scale  the  basis  functions  (2.1)  so  that 

max|K  (r.)|  =  1  Vs,j,  where  the  maximum  is  taken  over  all  collocation 
s  3 

points.  Denote  these  scaled  functions  by  f^(x,y),  f2(x,y) , . . . ,f_( x,y) ,  and 
3f . 

orthonormal ize  (x,y)  with  respect  to  the  inner-product 

m 

(f#g)  :=  I  f(x.,y.)  g(x.,y.),  (2.7) 

j=1  J  J  J  1 

£ 

where  we  sum  over  all  collocation  points,  and  —  denotes  the  outward  normal 

on 

/V 

derivative  w.r.t.  V.  We  denote  the  orthornormalized  functions  by  f^(x,y), 
f2(x,y>,...,f  (x,y)  .  We  seek  a  linear  combination 


*  ( x , y )  =  l  a  f  (x,y) 


(2.8) 


2 

and  constants  w  which  satisfy  (1.2)  and  (1.3)  at  the  collocation  points. 

This  gives  rise  to  an  overdetermined  generalized  eigenvalue  problem,  which  we 

T 

symbolically  write,  with  a  =  (a^,  a^,  •••>&  )  >  as 


m/2  eq.  from  top  nodes 


m/2  eq.  from  bottom  nodes 


(2.9) 


Multiply  (2.9)  by  the  transpose  of  the  orthonormal  left  hand  side  matrix.  The 
resulting  system  we  write  as 


2  2  2 
y  a  =  Aa  ,  u  =  1/v» 


(2.10) 


(2.10)  is  an  ordinary  eigenvalue  problem  with  an  n  x  n  matrix  A.  We  can 


regard  equation  (2.10)  as  Galerlein  equations  obtained  by  using  the  test 

3f1  3f2  3fn 

functions  T — ,  -r — ,  •••,  t —  and  inner-product  (2.7). 

3n  3n  3n 

2 

Solving  eigenproblem  (2.10)  yields  eigenvalues  y^  and  corresponding 

eigenvectors  a..  j=  (a.^..,  a2  j  '  ‘  *  ,an  j  ^  *  We  obtain  approximate  eigenvalues 
~2  2 

Wj  :=  1/y^,  to  problem  (1.1)-(1.5>  and  corresponding  approximate 


eigenfunctions 


<fi  . ( x,y )  :=  l  a  f  (x,y) . 
3  K=1 


(2.11) 


~2  ~2  ~2  ~2 

Enumerate  the  w.  so  that  w,  <  w,  <  ...  <  w  .  In  section  4,  we  derive  the 
]  l  i  n 

~2  2 

following  error  estimate  for  w^:  Let  w^  be  the  true  least  eigenvalue  for 
(1.1)- (1.5),  and  let  3Pfc  denote  the  part  of  3P  which  coincides  with  the  x- 
axis.  Let  3Pb  :®  3P  -  3Pfc.  Let  s  denote  arc  length.  Then 


~2  2  dV>  (x(s),y(s)) 

lwi  "  wil  *  J  - a -  +  ,(x(  s)  ,y(  s) ) )  |e  (s)|ds  +  (2.12) 

1  1  3pfc  as 

d<P  (x(s),y(s>)  2  2 

+  J - - -  +  (^  ( (x(  s)  ,y(  s)  )  )  |e  (s)|ds,  (2.13) 


V 


where 


t  ^ , 

efc(t)  :=  ]  (^j^-(x(s),y(s) )  -  w^  (x(  s )  ,y(  s ) ) )  (~j^(x(  s ) ,  y(  s )  ) )  ^s,  (2.14; 


c W> 


z  ar  av 

efa(t)  :=  J  (x(s),y(s))(— — (x( s ) ,y( s )  ) )  ds. 


(2.15) 


The  parameter  t  also  denotes  arc  lenth,  which  is  measured  from  the  vertex 
along  in  the  integrals  (2.12),  (2.14),  and  from  the  vertex  along  3D  ^ 
in  (2.13),  (2.15).  Especially,  t  =  0  corresponds  to  the  vertex. 

The  error  estimate  (2. 12)— (2. 15) ,  we  can  use  for  the  determination  of  a 


good  allocation  of  source  points  (x^,  y^).  In  all  computed  examples,  we 


obtained  good  accuracy  with  only  one  source  point  (x^,y^),  i.e.  in  (2.1)  we 


have  i  =  1,  and  we  restrict  the  description  to  this  case.  Generalizations 
to  l  distinct  source  points  are  immediate,  but  we  note  that  the  effort  to 
determine  them  increases  rapidly  with  i;  the  determination  of  l  source 
points  yields  a  2fc-dimensional  nonlinear  optimization  problem. 

Since  we  need  not  find  the  best  allocation  of  (x^/y^),  a  good  placement 
suffices,  we  can  use  a  simple  method  for  its  determination.  Choose  (x^y^) 
arbitrarily  but  fairly  near  the  vertex  at  D.  Determine  the  collocation 
points,  solve  the  eigenvalue  problem  (2.10),  and  evaluate  the  error  estimate 
(2. 12)  —  (2. 15).  Repeat  for  the  remaining  4  points  in  the  5-point-stencil,  with 
points 


{(x1#Y1 ), 


( x  1  +h ,  y  1 ) , 


(x1-h,y1 ), 


(x^y^h),  (x 


1'Y1 


-h)}. 


(2.  16) 


h  >  0  is  selected  so  that  the  whole  stencil  is  exterior  to  3V .  If  the 
error  estimate  is  smallest  at  one  leg  of  the  stencil,  say  at  (x1~h,y1),  then 
move  the  stencil  so  that  (x^~h,y^)  becomes  its  midpoint  and  compute  the 


error  estimate  at  stencil  points  where  not  known.  This  way  the  stencil  is 
moved  in  the  plane  until  the  smallest  error  estimate  is  obtained  at  its 
center.  In  this  case  the  stencil  is  shrunk,  i.e.  we  replace  h  by  h/2. 
Repeat  until  the  stencil  size  is  smaller  than  some  prescribed  6  >  0. 


3.  Numerical  Examples 


;<■- 


V.  V 


~2 


In  the  examples  we  denote  the  computed  eigenvalues  by  wR,  k=1,2,3,... 


~2  ~2  ~2 

ordered  so  that  <  ...  .  The  number  of  collocation  points,  we 


denote  by  m,  the  number  of  basis  functions  by  n.  For  wedges  with  straight 


m 


sides,  we  let  n  :=  m/2.  Increasing  the  ratio  —  does  not  change  the  order 


of  magnitude  of  the  error,  which  shows  that  the  distribution  method  for  the 


collocation  points  is  appropriate.  For  wedges  with  curved  sides  we  chose 


m 


~  *  3,  since  the  allocation  of  nodes  for  such  wedges  is  somewhat  arbitrary; 
n 


we  do  not  compute  the  conformal  map  from  the  exterior  of  the  unit  disk  to  the 
exterior  of  the  wedge,  but  use  instead  an  explicitly  known  map  from  the 


exterior  of  the  unit  disk  to  the  exterior  of  a  wedge  with  straight  sides. 


m 


Again,  increasing  the  ratio  ~  does  not  change  the  order  of  magnitude  of  the 


error . 

Let  a  denote  the  interior  angle  at  the  vertex  of  V .  All  computations  were 
initialized  by  choosing  the  source  point  (x^,y^)  =*  (-/2  cos(^),  /2  sin(~)), 
setting  m  and  n  to  the  smallest  values  in  the  tables  below  and  letting 
h  defining  the  stencil  (2.16)  be  1/2  for  example  3.1  and  1/4  for  the  other 
examples.  The  stencil  was  moved  as  described  in  the  end  of  section  2  until 
the  estimated  error  was  smallest  at  its  center  point.  This  point  we  below 
call  computed  source  point,  and  at  this  point  we  computed  the  eigenvalues  for 
larger  values  of  m  and  n  also. 

The  discretization  replaces  the  original  eigenvalue  problem  (1.1 )— (1.5) 


by  an  n-dimensional  one  that  has  n  eigenvalues,  which  are  approximated  by 


~2  ~2 

w,  =  w.  ,  j  =  1(1 )n.  All  eigenvalues  of  the  discrete  problem  need  not 
J  3*n 


correspond  to  eigenvalues  of  the  continuous  problem  (1.1)- (1.5),  and  we  need  a 


simple  criterion  to  determine  how  many  of  the  wj  n»  j  =  hDn,  to  accept  as 


approximations  of  eigenvalues  to  the  continuous  problem.  Our  criterion  is 


-11- 


•  /-’vV.VyWJ 


~2 

based  on  the  dependence  of  the  w  on  n.  For  some  integer  k,  depending 

J*n 

on  the  shape  of  "  but  independent  of  n,  the  w  ,  j  ■  1(1  )k,  converge 

3  »n 

as  n  increases,  as  far  as  this  can  be  judged  by  numerical  computations  with 

~2  ~2 

increasing  but  finite  n.  The  convergence  of  the  differences  w,  .  -  w , 

D-1#n  j  ,n 

~2 

as  n  ♦  ",  j  ■  2(1 )k,  is  as  fast  as  the  convergence  of  the  w.  ,  j  *  1(1 )k. 


For  j  >  k,  the  w.  decay  with  n,  but  limit  values  cannot  be  discerned, 

3/n 

~2  ~2  ~2  ~2 
only  that  w  >  w  ,  j  >  k.  The  pairwise  convergence  -  w .  to 

j  #n  k  rn  j#n 

0  as  n  ♦  j  >  k  is  very  pronounced,  however.  This  suggests  that  the 

~2 

k  smallest  computed  eigenvalues  w  ,  j  *  1 ( 1 )k ,  approximate  eigenvalues 

j  *n 

of  the  continuous  problem.  Numerically  we  determine  k  by  only  regarding 
~2  ~2 

differences  w,  .  -  w.  for  increasing  values  of  n.  We  choose  k  to  be 

j“i»n  j,n 


the  smallest  integer  such  that 


~2 


~2 

w,  «  -  w,  ♦  o  <  0,  n  ♦ 

k-1,n  k,n 


-2  ~2 

w .  „  -  w .  ♦  0,  n  ♦ 

j-1,n  j  ,n 


rj>k' 


where  5  is  any  constant  <0.  We  note,  in  passing,  that  the  ’eigenvalues’ 

~2 

w_j ,  j  >  k  are  more  sensitive  to  changes  in  computational  details,  like 

~2 

selection  of  source  point,  than  the  w  ,  1  <  j  <  k.  All  computations  were 
done  cm  a  VAX 11/780  in  double  precision  arithmetic,  i.e.  with  12  significant 
digits. 

Ex.  3.1.  We  illustrate  the  numerical  method  for  a  domain  with  known 
eigenvalues.  Let  D  be  a  wedge  with  straight  sides  and  an  interior  angle  of 
15  degrees.  Then  w^  =  sin  (~(2A-1>),  i=1,2,3,  are  the  eigenvalues  found  by 
Ursell  [9],  and  the  computations  also  identified  precisely  3  converging 

~2 

eigenvalues  by  the  criterion  described  above.  The  error  estimate  for  w1  is 


defined  by  (2. 12)— (2. 15).  The  table  holds  for  the  computed  source  point 


(x^y  )  -  (-6.351857,  2.30591). 


m 

n 

|w1-w1l 

error  est.  for  w^ 

Iw2-w2l 

|w3-w3l 

22 

11 

1.9»10~8 

8.4* 10”7 

i.o»io"4 

3.4* 10"2 

42 

21 

5.2*10_1° 

i.i* io”8 

2. 5» io”7 

ro 

1 

0 

• 

00 

e 

v— 

The  computations  suggest  that  there  are  no  other  eigenvalues  than  those  found 
by  Orsell  19].  e 

To  generate  various  bottom  shapes,  a  (p^,  p2»  p3 ) -family  of  curves  is 
used  and  for  convenience  defined  parametrically  in  the  complex  plane.  For 
t  <  0, 


”ip1 

x(t)  +  iy( t)  :»  (x^t)  +  iy^tjje  •  p2 

5t 

x^t)  s»  c»u  +  "^rd  ~  p3*u)  (3.1 

Yl(t)  :«  s*v  -  -f|{1  +  wv), 

where  c  :■  2  cos2(^),  s  sin(-~) ,  u  s*  tanh(-^), 
v  s*  tanh  ( 1  -  tan(^) ,  w  :■  (p3~l)tan(^)  -  1.  This  curve  forms  an 
angle  p^  with  the  positive  x-axis  at  the  origin,  x(0)  *  y(0)  =  0.  The 
parameter  p2  is  a  scaling  factor  introduced  to  make  t  approximately  equal 
to  arc  length  for  the  t-interval  under  consideration,  while  p3  determines 
the  shape  of  the  curve  before  it  asymptotes  a  straight  line  (as  t>-» ) . 

Ex.  3.2.  We  compute  eigenvalues  for  a  region  with  a  lower  boundary 
defined  by  a  curve  (3.1)  with  p,  -  tt,  p.  *  1,  and  p_  =  4,  see  figure 


3.1.  The  inscribed  (straight)  wedge  with  angle 
the  collocation  points. 


w 

12 


is  used  when  allocating 


Figure  3. 1 

~2  ^2 

k,  the  number  of  convergent  eigenvalues  w  is  2,  but  we  also  show  w,  as 

j  3,n 

being  a  typical  spurious  eigenvalue  of  the  continuous  problem.  The  computed 
source  point  is  (x^,y^)  ”  (-0.69500,  0.18459). 


m 

n 

-2 

W1 

~2 

error  est.  for  w^ 

-2 

W2 

-3 

W2 

-2  -2 
W3_W2 

32 

11 

0.305559 

—4 

9.4*10 

0.89902 

1.67 

0.77 

62 

21 

0.305434 

2.0*10“4 

0.86034 

1.26 

0.40 

92 

31 

0.305430 

3. 7* 10~5 

0.85668 

1.14 

0.28 

For  coeqparison,  we  note  that  the  eigenvalues  of  the  inscribed  (straight)  wedge 

—’2  f  T\ 

are  those  of  example  3.1.  w^  agrees  with  sin(— j  to  1  significant  digit.  ■ 
Ex.  3.3.  We  compute  eigenvalues  for  the  region  shown  on  figure  3.2.  Hie 
lower  boundary  is  defined  by  a  curve  (3.1)  with  p^  1  p2  ■  p^*  2.  The 
(straiqffit)  wedge  with  angle  ~  indicated  is  used  for  allocation  of 


-  ,*\'V  V  V  V  V  ■  V  v  V  v 


collocation  points. 


1 


The  computed  source  point  is  (x1fy^)  ■  (-1.74628,  -0.10797).  Two  converging 
-2  -2 

eigenvalues  w^,  w^  can  be  identified,  but  we  also  show  the  smallest  spurious 
~2 

eigenvalue  w^ . 


m 

n 

—2 

W1 

—2 

error  est.  for  w1 

** 

W2 

-2 

W3 

~2  -2 
W3_W2 

32 

11 

0.29001 

1.3*10-2 

0.8194 

1.58 

0.76 

62 

21 

0.28893 

2.0»10"3 

0.7746 

1.12 

0.34 

92 

31 

0.28888 

2.8*10"4 

0.7731 

1.06 

0.29 

/  IT  \ 

The  straight  wedge  on  figure  3.2  has  eigenvalue  sin[gj  *  0.342  and 

sin(~)  ”  0.866,  see  Ursell  [9],  After  rounding  to  1  correct  decimal 
/IT  \  -2 

sm^~J  and  w1  agree,  a 


4. 


Error  Estimate  by  Perturbation 


The  particular  problem  considered  can  be  embedded  in  a  slightly  wider 


class*  Suppose  V  is  a  domain  which  is  a  generalization  of  a  wedge  that  is 


semi-infinite,  with  two  smooth  (C;  )  boundaries,  a  'top'  (t)  and  a  'bottom' 


(b),  which  intersect  only  at  the  origin,  making  there  a  finite  non-zero  angle 


with  each  other.  These  boundaries  extend  to  infinity. 


The  slightly  more  general  eigenvalue  problem  is:  find  numbers  w  and 


functions  $  e  C  in  5  such  that 


A$  “  $  in  V 


(4.1) 


on  the  bottom  x  "  x^(s). 


(4.2) 


4^  m  w2$  on  the  top,  x  *  x  (s), 
an  "■  -"t 


(4.3) 


where  n(s)  is  the  (continuous)  unit  outward  normal  to  the  boundary,  and  s 


is  arc  length.  Also  ^  satisfies  a  finite  energy  condition  (1.5),  at  most 


♦(x)  ♦  o(|x|  *1+**),  6>0,  | x |  ♦  •  in  V. 


(4.4) 


For  these  problems  it  can  be  readily  shown,  by  using  standard  manipu¬ 


lations  together  with  Green's  Theorem,  that  the  eigenvalues  w  are  real  and 


2  2 


the  eigenfunctions  $  belonging  to  different  eigenvalues  w  ,  w  are 

3  X  j  K 


orthogonal#  in  that  both  the  boundary  arc  length  integral 


jt  ds  -  o,  j 


(4.5) 


and  the  area  integral 


.-  -  •'v\  YvV\  . 


.  v-;.- ' 

\  V  --  ■  •-  •- 


\  \  \  ^  v  %  . 


jyv* +  *  -V  **  ■  °'  3  * k 


(4.6) 


vanish.  Because  the  relation 


w 2  -  I] 
j  JJV 


[IV^I2  +  ♦*!  as 
jt*2f 


(4.7) 


between  eigenpairs  w^,  also  holds,  it  follows  that  the  eigenvalues  are 

2 

all  positive,  so  the  notation  w^  is  justified. 

An  alternative  foraulation  of  the  problem  for  the  least  eigenvalue  now 
follows.  Minimise  the  functional 


"  JipWI2  +  *2]ds 


2 

over  all  ♦  €  C  in  tJ,  which  satisfy  the  normalizing  relation 


(4.8) 


Jt*  da  -  1, 


(4.9) 


which  vanish  in  $  sufficiently  rapidly  that  the  two  integrals  (4.8,  4.9)  are 


defined. 


The  calculus  of  variations  then  reveals  that  a  minimizing  function  4> 


satisfies 


At  *  ♦ 


0 

dn 


in  V 


on  the  bottom 


on  the  top 


(4.10) 


where  X  is  the  Lagrange  multiplier.  By  virtue  of  the  preceeding  discussion 


X  is  also  the  value  of  Ifl  at  the  minimum,  and  hence  the  least  eigenvalue. 


A  Perturbation  Analysis 

The  formulation  (4.8-4.10)  above  allows  some  qualitative  statements  to  be 
made  for  problems  in  which  the  boundaries  of  V  undergo  perturbations.  The 
aim  is  to  estimate  the  change  in  the  least  eigenvalue  for  a  prescribed 
perturbation . 

Let  the  top  boundary  perturbation  t*  be  described  by  an  equation 


x(s)  *  x  (s)  -  e  (s)n(s) 
-  t  - 


et(0)  -  0 


(4.11) 

(4.12) 


where  the  boundary  perturbation  function  e  is  everywhere  positive  or  zero. 


and  together  with  its  derivative,  is  small 

de. 


le  I,  It-**  <  E  <  1, 
t  ds 


where 


I f I  -  sup  | f ( s ) | . 

0<s 

Denote  as  p*  the  region  interior  to  the  perturbed  top  and  the  original 


bottom. 


Suppose  w2  is  the  least  eigenvalue  for  the  problem  (4. 1-4.4)  on  Q 


with  corresponding  eigenfunction  $,  and  let  0  be  the  least  eigenvalue  for 


2  2 

the  problem  on  P* .  Since  V  C  t? ,  one  would  anticipate  that  0  <  w  ;  the 


proof  that  this  is  indeed  the  case  follows. 


Use  as  a  trial  function  to  estimate  0  the  original  eigenfunction  +, 
which  must  first,  however,  be  renormalised.  Thus  one  has  the  formula 
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,  *\ 
V. 


«r.  *  -  ■ 
V  v/ 


(4.13) 


"  *  “  y*‘£t  d8t  +  °<e2>) 


/t(*2-2£tw2*2)dst(l  +  0(6^1) 
(l-2w2  Jt  et*2ds)(l  +  o(e2)) 


after  using  tha  first  term  in  a  local  Taylor  expansion  of  +  on  t  to 
estiMte  its  value  on  t',  the  boundary  condition  (4.4),  and  after  observing 
that  the  arc  length  eleMnts  on  t*  and  t  satisfy 

dsfcl  *  dsfc(l  +  0(e2)) 
by  the  aisused  properties  of  c^. 

So  the  variational  result  (4.8)  yields 


lnt(*> 


The  quantity  Ip_p,(4)  can  be  estimated  in  the  limit  6+0 


Ip_p<*>  “  /t(|V*12  +  *2)  6t  68  +  °(e2) 


(4.14) 


(4.15) 


/t  ((“^)2  +  (Hw4)*2)  ds  +  0(e2 


so  that  on  combining  the  three  formulae  (4.13,  4.14,  4.15)  it  follows  that 


G2  <  w2(l-w  2jt<^)2  etds-(w2+w  2>jt$2et<Ss  +  0(e2  >)(  1-2w2J ^^ds)”1  ( 1+0(e2 )) 


(4.16) 


w2(l-v  2Jt^)2  etds  "  (w  2-w2)Jt$2et<is)  ( 140(C2  ) ). 


Provided  0  <  w  <1  (as  is  the  case  in  the  computed  examples)  the  factor 
-2  2 

(w  -  w  )  is  positive,  and  because  >  0  by  hypothesis,  it  follows  that 


n2  2 
ft  <  w  . 


(4.17) 


Por  a  bottom  perturbation  which  yields  a  perturbed  region  V '  C  V  let 


x(s)  “  x.  (8)  -  e. (s)n 


(4.18) 


define  the  perturbed  bottom  b*,  where  e  has  properties  similar  to  those  of 

o 


e^,  that  is 


le.l ,  Irr8!  <  e  0,  e.(0)  -  0 

D  08  DO 


(4.19) 


As  a  trial  function  4  to  estimate  the  perturbed  least  eigenvalue  ft 
of  V  the  corresponding  eigenfunction  4  of  P  is  used.  Then  the 
variational  formulation  requires  that 


ft  <  V  ♦ 


'  VP****) 


(4.20) 


Hence  the  result 


2  2 
ft  4  W 


(4.21) 


again  follows. 


Thus,  in  summary,  when  P*  is  obtained  from  P  by  either  (or  both) 

boundary  perturbation(s)  of  the  type(s)  described,  (and  hence  P'  c  V),  the 

2  2 

least  eigenvalues  fl  ,  w  of  P'  and  P,  respectively,  satisfy 

2  2 

fl  <  w  .  (4.22) 

The  converse  is  also  true,  since  the  boundary  perturbation  process  described 
is  reflexive.  That  is,  if  P'  2  P»  then 

fl2  >  w2.  (4.23 ) 

While  the  results  (4.22,  4.23)  can  be  used  for  a  priori  estimation  using 
known  solutions  as  comparison,  the  relationship  between  the  least  eigenvalues 
of  P,  and  P'  (obtained  by  perturbation)  can  be  refined.  First  suppose 

2 

P  *  2.  P  and  the  eigenfunction  4  corresponding  to  the  least  eigenvalue  w 
of  P  is  known.  Extend  the  eigenfunction  4  to  P*  by  constructing  a 
preliminary  extension  4^  which  coincides  with  4  in  P  and  has  4.| 
constant  along  any  given  normal  n'  to  t',  b* ,  taking  the  value  on  that 
normal  equal  to  the  value  of  4  where  the  normal  intersects  t.  (Multi¬ 
valuedness  is  excluded  by  the  smoothness  assumptions  on  e^,  e^).  Then  4>  1 

2 

can  be  modified  to  produce  a  suitable  trial  function  4  e  C  in  P '  by 

2 

introducing  large  curvature  (“—r)  in  suitably  small  transverse  neighborhoods 

dn 

of  t,  b.  Such  a  process  can  be  carried  out  so  that,  for  any  desired  m 

Ip(4)  -  Ip(4>(l  +  0(em) ),  (4.25a) 

and 

ds  *  ds^  +  *  “  4(^  +  0(em)) ,  on  t,  and  b.  (4.25b) 

Having  constructed  this  4,  one  notes  that  its  normalization  is  essentially 
unchanged  since 
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it,*2dst'  =  It*2*8^1  +  0(e;2)) 


(4.26) 


by  the  properties  of  the  perturbation. 


Since  V'  ^  V  the  inequality  (4.23)  requires  that 


2  „  2 
w  <  ft 


and  the  minimum  principle  (4.8)  then  demands  that 


w 2  <  n2  <  Ip.(*> 


or 


(4.27) 


0  <  n2  -  W2  <  ip(*)  +  ip,_p(*)  -  v2. 


But  equation  (4.25)  and  the  definition  of  ♦  show  that 


2  m 

Ip(*)  -  w  =  0(e) 


(4.28) 


so  the  basic  inequality  (4.27)  becomes  essentially 


0  <  n2  -  w2  <  !«,«(*). 


V'-V 

In  the  limit  as  e+0 ,  the  right  hand  side  of  the  last  inequality  can  be 


(4.29) 


computed  by 


1V'-Vli)  =  ^t^17*'2  +  *2)|et|dst  +  ib^17*'2  +  *2^|Eb|dsb  +  °(eZ)  (4,30) 


d« 


where  V4,t  are  evaluated  in  the  interior  of  (P '  -P).  Here  “  vanishes 

an 


by  construction,  or  =  0, 


!„.„(♦>  *  Ut  +  *  ♦2)l£*lds.)  ('  *  0<e2>) 

-  (it  +  Jb(©2  *  (1  £0(c2)) 


(4.31 ) 


where  the  suffix  *  denotes  t  or  b  as  appropriate. 

Thus,  in  the  situation  described  the  results  (4.27),  (4.31)  can  be 

summarized  as 
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0  <  f22  -  w2  <  Ip, _p( ♦ )  (4.32) 

and  I p, _p ( ♦ )  can  be  approximated  by  the  {  }  factor  in  equation  (4.31). 

One  can  proceed  analagously  to  estimate  the  relation  between  the  least 

2 

eigenvalues  of  P,  P*  when  P'  ^  P  and  the  eigenpair  w  ,  $  belonging 
to  P  are  known.  It  will  follow  that  a  result  similar  to  the  above  (4.32) 


holds,  namely 


0  <  w2  -  f)2  <  1(f) 


(4.33) 


after  a  calculation  which  is  similar  to  the  one  used  in  developing  the 
inequality  (4.32),  but  different  in  detail. 

When  p'  is  neither  wholly  interior  to  P,  nor  contains  it,  a  two  step 
procedure  is  carried  out  to  relate  the  least  eigenvalue  w  of  P  to  that 
of  P'.  Denote  as  P*  the  region  P  U  (P'-P),  and  use  as  a  trial  function 
in  both  V  *  and  P'  the  eigenfunction  $  of  P  extended  by  the  rules 
previously  described:  that  is,  where  t',  b'  fall  outside  P,  extend  $  as 
a  constant  along  normal  lines  to  these  boundaries,  taking  a  value  on  any  such 
line  equal  to  the  value  of  $  where  that  line  strikes  t,  or  b.  At 
boundaries  t,  b,  $  is  adjusted  to  satisfy  smoothness  conditions  or  boundary 
conditions  as  necessary,  but  in  regions  whose  transverse  extent  is  arbitrarily 
small  in  comparison  with  e.  Noting  that  P*  D  P»  P*  3  P1,  one  deduces  the 


pair  of  inequalities 


0  <  a2  -  n2  <  I  *  ,(♦) 

V  -V 


0  <  A2  -  w2  =  I  #  (<|>) 

V  -v 


(4.34) 


from  the  properties  of  the  trial  function,  where  A  is  the  least  eigenvalue 
of  P*.  From  these  inequalities  (4.34)  it  follows  that  a  conservative 
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estimate  is 


o  <  |w2  -  n2|  <  +  ip*_pI<*> 

and  the  last  two  integrals  can  be  approximated  as 


(4.35) 


W*’ 


W<4‘ 


/tt( 


ail 

ds' 


+  f  | 


£tlld8 


JbHdsJ 


+  <f>  |e  lids  (4.36) 


The  inequality  (4.35)  and  the  evaluation  (4.36)  form  the  basis  of  the 
computable  error  estimate. 


An  Error  Estimate 

In  estimating  the  error  in  the  computation  of  the  least  eigenvalue,  the 
premise  is  that  the  small  residuals  in 

to-V*'  <4-37’ 

and 

i?i  _  «  d  (s)  (4.38) 

dn  t 

2 

computed  on  b,  t  respectively,  the  boundaries  of  0  suggest  X  ,  $  are  an 

eigenpair  of  the  least  eigenvalue  of  a  region  O' t  which  is  a  perturbation 

of  0  •  (By  construction,  $  is  a  solution  of  equation  (4.1)  in  a  region 

including  0) •  Knowledge  of  the  $  field,  and  the  location  of  b,  t  enable 

the  relative  location  of  b',  t'  to  the  established,  and  hence  (using  the 

.22  2 

results  of  the  perturbation  analysis)  a  bound  on  | y  -  X  |  where  y  is  the 
true  least  eigenvalue  of  p. 

In  the  case  of  the  bottom  boundary,  for  example,  for  the  given  $  field 

it  is  necessary  to  determine  the  nearby  curve  b' 

x  =  x.  (s)  +  e  (s)n(s)  (4.39) 

~  ~b  b  ~ 

on  which 

•  n'  =  0  (4.40) 
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where  jj*  *  n  is  the  local  unit  outward  normal  to  b'.  When  the  curvature 
of  b  is  not  large,  s  and  n  coordinates  may  be  regarded  as  Cartesian,  the 
boundary  condition  becomes  approximately 


d&  de,  d$ 

fo  ^  ^  0 

ds  ds  dn 

on  b',  or,  approximately  (equation  (4.37)) 


Thus  the  ultimate  expression  is 

Vs)  =  -Jo  V§ 


ds. 


x(  s) 


(4.41  ) 


(4.42) 


(4.43) 


allowing  for  the  boundary  condition  £^(0)  =  0.  Similarly,  on  the  top 
boundary  one  requires 

0  -  V$.n*  -  y2* 


(4.44) 

*  -srS*1  * 0(c>)  * 


2  2 

(assuming  the  calculations  are  good  enough  that  y  is  close  to  X  ),  so 
that  taking  the  local  approximate  values 


it  follows  that 


d$ _ 

dn' 


d£ 

ds 


x+e.  n 


(4.45) 


*.  “  J’  <vl*l 


ds 


(4.46) 


ppigrBTC  ■■  a  -v  t v*  mw.-.-'.-v 


,Sfr 

>-£ 


>.v 
• .  •  * 

r>s- 

s.  ,\  . 
■.■„•••. 


v.'J 


The  results  {4.43,  4.46)  can  then  be  used  In  the  formulae  (4.35,  4.36)  to 
estimate  the  error  in  the  computation.  (The  formulae  (4.43)  and  (4.46)  are 
constructed  in  the  knowledge  that  ~^L  .  vanishes  only  in  the  limit  I x| 

Q8  t,D  ~ 

when  t  i b  the  eigenfunction  corresponding  to  the  least  eigenvalue.) 
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